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Abstract 

We use the T-matrix approximation to analyze the effect of a localized impurity on 
the local density of states in mono- and bilayer graphene. For monolayer graphene the 
Friedel oscillations generated by intranodal scattering obey an inverse-square law, while 
the internodal ones obey an inverse law. In the Fourier transform this translates into 
a filled circle of high intensity in the center of the Brillouin zone, and empty circular 
contours around its corners. For bilayer graphene both types of oscillations obey an 
inverse law. 

Graphene has been studied extensively in the recent years. Its most fascinating aspect is 
the existence of linearly-dispersing gapless excitations in the vicinity of the Dirac points. This 
gives rise to very interesting electronic properties such as Friedel oscillations in the local density 
of states at low energy which decay as 1/r 2 [HE], instead the usual 1/r characteristic to two- 
dimensional systems (3j H] . It is very important to understand well the physics of impurity 
scattering in graphene, by studying for example the local density of states (LDOS) in the 
presence of single impurity scattering [TJ [2j [3], [5] . Comparison with experiments can provide 
information about graphene's fundamental physics, and about the nature of the impurities. 

In this Letter we analyze the Fourier transform of the density of states measurable by 
Fourier transform scanning tunneling spectroscopy (FTSTS). Such measurements have re- 
cently been developed for graphene [6], as well as for other two-dimensional materials such as 
ErS%2 [7] and the cuprates [8]. 

Our first observation is that the FTSTS spectra can be used to distinguish between mono- 
layer and bilayer graphene. In particular, for monolayer graphene with a localized (delta- 
function) impurity potential, the Friedel oscillations generated by intranodal scattering decay 
as 1/r 2 at low STM bias, consistent with previous analysis [HE]- In the FTSTS spectra this 
is manifested by a filled circle of high intensity in the center of the Brillouin zone (BZ), with 
a radius proportional to the STM bias. On the other hand, the Friedel oscillations generated 
by the scattering of quasiparticles between different Dirac points decay as 1/r. In the FT- 
STS spectra these oscillations are translated into circular contours of high intensity centered 
around the corners of the BZ and around sites of the reciprocal lattice. Due to the form 
of the underlying Hamiltonian, the distribution of intensity on some of these circles is not 
rotationally invariant. 



For the bilayer system, at low energy the oscillations due to both intranodal and internodal 
scattering have a 1/r dependence, corresponding to circular lines of high intensity close to 
center and the corners of the BZ, and around the sites of the reciprocal lattice. At higher 
energies the splitting of the bands for the bilayer sample is also observable in the FTSTS 
spectra. 

We also note that the FTSTS spectra can distinguish between different types of impurities. 
For example, for the case of a screened-charge impurity, the effect of internodal scattering is 
greatly reduced compared to the effect of intranodal scattering. This gives a clear signature 
in the FTSTS spectra which can be observed in an experiment. 

Our last observation is that the FTSTS spectra, besides providing information about 
the band structure of graphene, can also give insight into the underlying Hamiltonian. In 
particular, the shift of the decay of the Friedel oscillations from 1/r to 1/r 2 , and the rotational 
asymmetry of some of the high-intensity spots, are strongly dependent on the peculiar form of 
the tight-binding Hamiltonian, and cannot be deduced solely from band-structure arguments. 

The tight-binding Hamiltonian for monolayer graphene is: 

n = J d 2 k[aib % f{k) + h.c], (i) 

where the operators a\ correspond to creating electrons on the sublattice A and B re- 
spectively, and f{k) = —t Y^j=i exp(ifc • a,-). Here a\ = a(y/3x + y)/2, a 2 = a(— V3x + y)/2, 
a 3 = —ay, t is the nearest-neighbor hopping amplitude, and a is the spacing between two 
adjacent carbon atoms, which we are setting to 1. 

We will use this form of the Hamiltonian when performing our numerical analysis of the 
FTSTS spectra. However, it is useful to expand the Hamiltonian close to the corners of the 
BZ, which we also denote as nodes or "Dirac points," and use the linearized form to solve the 
problem analytically at low energies. The momenta of the six corners of the Brillouin zone 
are given by K h2 = [±4^/(3^,0], K 3>4 = [±2tt/(3V3), 2tt/3], K 5>6 = [±2tt/(3>/3), -2tt/3]. 
Close to each corner, m, of the BZ we can write f(q + K m ) rs (f) m (q) = v m q ■ J m , where q 
denotes the distance from the respective corner. Also v\p = St/2 = v, = v exp(— in/S), 
r 5j6 = v exp(«7r/3) and J 1)2 = (±l,-i), J 3A = J 5;6 = (±l,i). 

The corresponding Green's function, Q(k,u), derived from the tight-binding Hamiltonian 
in Eq. ([1]) can be expanded at low energy around the six nodes (denoted m), and in the 2x2 
(A,B) sublattice basis can be written as: 

g(k^)^G m (k^) = L^( u + * Mk)\ (2) 

where 8 is the quasiparticle inverse lifetime. The Fourier transform of the linearized Green's 
function is given by: 

n <- \ ( H o 1] ( z ) ^0 m (r)^ 1 (1) (^)^ ( ^ 

where z = uor/v, H^\(r) are Hankel functions, r = |f|, and <f) m (r) = v m r '■ J m /(vr). 

We first focus on a delta-function impurity localized on an atom belonging to sublattice 
A. In the (A, B) basis the impurity potential matrix V has only one non-zero component 
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Vii = u. We start with a T- matrix study [TOl [3] of the full Hamiltonian presented in Eq. (flT) , 
and we analyze our results numerically for various energies. The resulting FTSTS spectra 
(corresponding to the real part of the Fourier transform of the LDOS) are plotted in Fig. [TJ 
There are several interesting features that should be noted. First, there are regions of high 
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Figure 1: FTSTS spectra for a monolayer graphene sample with a single delta- function impurity. 
Figs, l.a), l.c) and l.d) correspond to energies 0.15i, 0.6i, and 1.2t respectively, at <5 = 0.07i. The 
BZ is indicated by dashed lines. The actual lowest (0) and highest (1) values of the FTSTS intensity 
are different for each energy ((—1.3,2.6) for 0.15i, (—0.8,5.9) for 0.6i, and (—6.2,7.2) for 1.2i in 
arbitrary units). Fig. l.b) shows a cross section of the FTSTS intensity as a function of k x for 
k y = 0, and for energy 0.15t. 

intensity in the FTSTS spectra corresponding to intranodal quasiparticle scattering (central 
region) and internodal scattering (outer regions). Similar features have also been observed 
experimentally [6]. The high intensity regions that we find are point-like at zero energy, and 
acquire distinct features as one increases the energy (STM bias). Thus, at low energy the 
central high-intensity region is a filled circle, while the outer regions are empty. Also, the 
rotational symmetry of the high-intensity regions located at the corners of the BZ is broken, 
while it is preserved for the high-intensity regions centered on sites of the reciprocal-lattice. 
With increasing the energy even higher, other effects such as the changing of the shape of the 
equal-energy contours from circular to triangular (trigonal warping) start playing an important 
role. At very high energy the FTSTS intensity map becomes quite intricate. 

We now turn to the analytical study of the dependence of the LDOS on the relative position 
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with respect to the impurity (r) at low energies. In this range, the physics is dominated by 
linearly dispersing quasiparticles close to the Dirac points. We find that spatial variations of 
the LDOS due to the impurity are given by: 

p(f, E) oc -Im[£(-f, E)T(E)Q(r, Im[e^-^G m (-f, E)T(E)G n (f, E)}, (4) 

m,n 

where m,n denote the corresponding Dirac points. Here T(E) is the T-matrix, which for a 
delta-function impurity is given by [10J T(u) = [I — V J ^^Q(k,uo)]~ l V, where / is the 2x2 

identity matrix, and the integral over k is performed on the BZ, whose area is Sbz = 87r 2 /3\/3. 

Using Eq. ([3D and expanding the Hankel functions to leading order in 1/r, we find that far 
from the impurity (uir/v ^> 1) the corrections to the local density of states due to scattering 
between the nodes m and n are given by: 

p mn (r, u) oc ^Imitiujy^-^+^i [l - C(r)0 n (rj] }■ (5) 

where t(u) is the non-zero element of the T-matrix (Tn), and we used the fact that 4>{—r) = 

-0(f). 

In the case of intranodal scattering (m = n) the above expression vanishes and the LDOS 
is dominated by the next leading correction p m (f,u) oc sin(2u;r/i> )/r 2 . This is different from 
what usual wisdom would suggest for a two-dimensional system (1/r decay) [31 H], and has also 
been described in Refs. [HE]. We should note that the two-dimensional FT of sm(2ur/v)/r 2 
is roughly p m (q,u>) oc tt9(2uj — qv)/2 + arcsin(2o;/gw)[l — 6{2u — qv)]. This corresponds to a 
filled circle of high intensity in the FTSTS spectrum, which is consistent with the results of 
our numerical analysis for the central region of high intensity. 

Nevertheless, for the decay of the Fridel oscillations generated by internodal scattering 
(m ^ n), the leading order behavior is 1/r. The FT of cos(2ujr/v)/r is 6(qv—2co) / \J q 2 v 2 — 4uj 2 , 
which translates into empty circles of high intensity in the FTSTS spectra, consistent with 
our numerical analysis. However, since the inverse quasiparticle lifetime S is finite, there will 
be some broadening of the resonances and some weight inside the circular contours. 

For some of the Fridel oscillations generated by internodal scattering, such as the ones 
between neighboring node pairs (e.g. (1,3) ), the rotational symmetry is broken: 
(fil(r)(j) 3 (r) = e~ ln / 3 (x + iy) 2 . However, for next-to-nearest-neighbor node pairs (e.g. (3,5)), 
03(r)0 5 (r) = e~* 2?r / 3 , and the oscillations are rotationally invariant. This is consistent with the 
results of our numerical analysis: the high-intensity regions centered on sites of the reciprocal- 
lattice are rotationally symmetric, while the high-intensity regions close to the corners of the 
BZ are not. 

We now switch gears and consider the case of a single impurity in bilayer graphene. The 
bilayer graphene consists of two graphene layers stacked on top of each other such that the 
atoms in the sublattice A of the first layer occur naturally directly on top of the atoms in the 
sublattice B of the second layer [11], with a tunneling coupling of t p . We consider the case 
of an impurity located on the sublattice A. The case of a single impurity located on a site of 
a different type, as well as the case of multiple impurities will be presented elsewhere. The 
resulting FTSTS spectra for the LDOS in the top layer are presented in Fig. 

Note that there are similarities and discrepancies between the monolayer and bilayer cases. 
Like in the monolayer case, there are areas of high intensity centered on the corners of the 
BZ, as well as on the sites of the reciprocal lattice. The main difference at low energy is 
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Figure 2: FTSTS spectra for a bilayer sample. Fig. 2. a) and 2.c) depict the FTSTS intensity in 
arbitrary units at energies O.li, OAt, with t p = 0.3t, and 5 = 0.05t. The actual lowest (0) and highest 
(1) values of the FTSTS intensity are (-6.9,3.6) and (-6.0,7.4) respectively. Figures 2.b) and 2.d) 
are cross-sections at k y = of Figs. 2. a) and 2.c). 
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that the central region of high intensity is an empty circle, and not a full circle (as for the 
monolayer case). At high energy, we also note a doubling of the number of high intensity lines 
corresponding to the doubling of the number of bands. 

An analytical study can be performed at low energies starting from the expansion of the 
Hamiltonian around the Dirac points m. In the sublattice basis (A,B) this yields [12] : 

n^(k) = [ fc ^ )]2 [ ^ )]2 | (6) 

where for simplicity we have set the effective mass of the quadratic spectrum to 1. The 
corresponding Green's function in real space is given by: 

G ^*\-u, {f)fH?(z) H m (z) < 7 > 



where we denoted z = ry\u)\/v. Starting from Eq. (|4j), we perform a similar analysis to the 
case of monolayer graphene. Thus we note that at large distances (z> 1), as opposed to the 
monolayer case, the leading (1/r) contribution for intranodal scattering is non- vanishing: 

p m (r,u) oc — \=cos(r>f\tj\/v). (8) 
\u\ 



This is consistent with the appearance of an empty circular contour at the center of the BZ, 
as opposed to the filled circle for the monolayer case. The leading contribution to the decay 
of the oscillations due to internodal scattering is also 1/r. Note that in the monolayer case 
the amplitude of the Friedel oscillations corresponding to intranodal scattering is independent 
of energy, and the amplitude of the Friedel oscillations corresponding to internodal scattering 
increases as uj 2 , while for the bilayer case all oscillations decrease with energy as l/y/\u\. 

The third system we consider is a monolayer graphene sample with a single screened- 
Coulomb (charged) impurity. Its impurity potential (in momentum space) has V\\ = V22 oc 
l/(9+l/eo), where eo ~ 4.5a is the screening length [13]. Given the momentum dependence of 
the impurity, one can no longer use the T-matrix approximation, but the Born approximation. 
Our results are plotted in Fig. Note that the intensity of the outer areas of high intensity 
is very much reduced compared to the intensity of the inner area. This is consistent with the 
form of the scattering potential, which generates less scattering between quasiparticles located 
on different nodes than between quasiparticles located on the same node. Similar qualitative 
differences have also been found between the LDOS oscillations generated by the "mirrage" 
images of Coulomb-type impurities and delta-function impurities in graphene PN junctions 

P-0- 

To conclude, we have computed the effect of single-impurity scattering on the LDOS, in 
particular on the Fourier transform of the LDOS measurable experimentally by FTSTS. We 
have found that the FTSTS spectra in the vicinity of an impurity are a very good tool to dis- 
tinguish between monolayer and bilayer graphene. In particular, for monolayer graphene, the 
Friedel oscillations due to intranodal scattering decay as 1/r 2 and are rotationally invariant. 
In the FTSTS spectra they correspond to a filled high-intensity circular region at the center of 
the BZ. On the other hand, the Fridel oscillations due to internodal scattering decay as 1/r. 
They lead to different FTSTS features, such as circular contours of high intensity centered 
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Figure 3: FTSTS spectra for a graphene monolayer with a screened Coulomb impurity at energy 
0.2t, and 5 = 0.07t. The lowest and highest values of the intensity are (—1.7, 42.3) in arbitrary units. 
Fig. 3.b) presents a cross-section at k y = 0. 



around the corners of the BZ and on sites of the reciprocal lattice. Some of these contours 
display a breaking of rotational invariance. For the bilayer case both the internodal and the 
intranodal Fridel oscillations decay as 1/r. 

We have also showed that the FTSTS spectra can be used to distinguish between different 
types of impurities, for example between a delta-function impurity and a screened Coulomb 
scatterer. 

Last, but not least, we have noticed that, while the FTSTS spectra in the presence of 
an impurity can give information on the band structure, they are not fully determined by it, 
but also contain very important information about the specific form of the Hamiltonian. We 
believe that this feature is very important, and could be also used in the case of cuprates 
to understand the physics of high temperature superconductivity. It would be interesting for 
example to compare ARPES measurements and FTSTS spectra in the presence of an impurity. 
For graphene, this would clearly distinguish the effects of the Hamiltonian from the effects of 
the band structure. For other systems such as the cuprates the differences between ARPES 
and FTSTS will be more complex, and may provide even deeper insights into the underlying 
physics. 
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